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Abstract —The development of hydraulic fracturing technology 
has dramatically Increased the supply and lowered the cost 
of natural gas in the United States, driving an expansion of 
natural gas-fired generation capacity in several electrical inter¬ 
connections. Gas-fired generators have the capability to ramp 
quickly and are often utilized by grid operators to balance 
intermittency caused by wind generation. The time-varying 
output of these generators results in time-varying natural gas 
consumption rates that impact the pressure and line-pack of the 
gas network. As gas system operators assume nearly constant 
gas consumption when estimating pipeline transfer capacity and 
for planning operations, such fluctuations are a source of risk 
to their system. Here, we develop a new method to assess this 
risk. We consider a model of gas networks with consumption 
modeled through two components: forecasted consumption and 
small spatio-temporarily varying consumption due to the gas- 
fired generators being used to balance wind. While the forecasted 
consumption is globally balanced over longer time scales, the 
fluctuating consumption causes pressure fluctuations in the gas 
system to grow diffusively in time with a diffusion rate sensitive 
to the steady but spatially-inhomogeneous forecasted distribution 
of mass flow. To motivate our approach, we analyze the effect 
of fluctuating gas consumption on a model of the Transco gas 
pipeline that extends from the Gulf of Mexico to the Northeast 
of the United States. 

Index Terms —Natural Gas Networks; Gas-Electric Coupling; 
Stochasticity; Reliability 

1. Introduction 

A dominant new load on gas pipeline systems is natural 
gas-fired generators m, El. An example of this dramatic 
change is seen on the gas pipelines that supply the elec¬ 
trical grid controlled by the Independent System Operator 
of New England (ISO-NE) where natural gas-fired electrical 
generation increased from 5% of total capacity to 50% in a 
span of 20 years 0. A parallel development in many U.S. 
electrical grids is the expansion of intermittent renewable 
generation such as wind and photovoltaic (PV) generation—a 
trend that is expected to continue as utilities work to meet 
renewable portfolio standards a, Q that mandate a certain 
fraction of electrical generation be derived from renewable 
sources. In contrast to traditional coal, hydro or gas-fired 
generation, these intermittent renewable generators have lim¬ 
ited controllability. To maintain balance of generation and 
load, other grid resources must respond to counteract these 
new fluctuations. Although many different types of advanced 


control of nontraditional resources are under consideration to 
provide balancing services, e.g. grid-scale battery storage and 
demand response, the control of fast-responding traditional 
generation (i.e. gas) is the current state-of-practice. 

Gas pipelines have traditionally supplied Load Distribution 
Companies (LDC) that primarily serve space or water heating 
loads that evolve slowly throughout the day in a relatively 
well-known pattern that is predicted based on historical in¬ 
formation and weather forecasts. Other traditional pipeline 
customers are industrial loads that change from day-to-day, 
but are very predictable over the span of twenty-four hours. 
The combination of expanded natural gas-fired generation and 
its use to balance intermittent renewable generation is creating 
loads on natural gas pipelines that are significantly different 
than historical behavior and will challenge the current pipeline 
operating paradigm that is used to control gas pressure. 

The flow in natural gas pipeline is determined via bilateral 
transactions between buyers and sellers in a day-ahead market 
with market clearing and gas flow scheduling done in advance 
of the subsequent 24-hour period of gas delivery. Scheduling 
consists of determining the locations and constant rates of gas 
injections. The initial market clearing assumes that gas con¬ 
sumptions are uniform over the subsequent 24-hour delivery 
period. Over the gas day, gas buyers improve their estimate 
of actual gas needs, and mid-course corrections are allowed 
through the transaction and scheduling of gas flows in two 
subsequent intra-day markets at 10 and 14 hours after the start 
of the 24-hour delivery period. 

When serving traditional gas loads, the variability during the 
gas day is relatively small and slow and is well managed 
by linepack, i.e. compressed gas stored in the pipeline. The 
pressure in a gas transmission pipeline ranges between a 
maximum set by engineering limits and a minimum delivery 
pressure set by contracts. A typical maximum pressure is 
around 800 psi, and flow of the gas causes the pressure to 
fall along the pipeline. As the minimum pressure 500 psi) 
is approached, gas compressors installed along the pipeline are 
used to boost the pressure back near the maximum. Typical 
spacing between compressors is ~ 50-100 km 
The relatively high operating pressures enable large gas trans¬ 
fer rates, and the spread between maximum and minimum 
pressure allows the pipeline to operate with an imbalance 


of gas injections and consumptions for hours at a time. An 
injection-consumption imbalance modifies the amount of gas 
stored in the pipeline via pressure changes, i.e. changes to 
the linepack. Linepack is sufficient to buffer the imbalance 
when serving traditional gas loads. However, the hydraulic 
fracturing-driven expansion of natural gas-fired generation 
capacity and its use to balance intermittent renewable 
fluctuations will result in larger and faster fluctuations in 
consumption (and possibly production) creating challenges to 
historical pipeline operations and reliability. 

The analysis in this manuscript is motivated by these new 
challenges. Our approach is built on top of any solution 
of the steady gas flow problem that determines the spatial 
dependence of gas flow and pressure and the dispatch of gas 
compressors to maintain pressure. For example, the steady 
flow solution can be found by solving an optimal gas flow 
(OGF) problem Cl or using a model that approximates com¬ 
pressor dispatch decisions in current gas pipeline operations. 
Using these steady solutions, we build on the ideas of m and 
develop analysis tools to provide a probabilistic measure of the 
impact on pipeline reliability created by stochastic deviations 
of gas consumption from the forecasted values used in the 
steady solution (and scheduled during market clearing). These 
new tools are based on a linearization of the basic gas flow 
equations around the forecasted solution that retains the effect 
of stochasticity in consumption. The effect of this stochasticity 
is assessed on a model of the Transco gas pipeline that extends 
from the Gulf of Mexico to the Northeast of the United States 
(See Fig. and Cl)- 

This manuscript builds on recent work 0, la to develop a 
theoretic and computational approach to analyze the evolution 
of pressure in a gas system over time and space when the 
system is imbalanced. The three main contributions of this 
approach are: 

• An analysis of spatiotemporal behavior of line-pack when 
the pipeline is subjected to stochastic gas consumptions. 
We observe that, even when fluctuations of the consump¬ 
tion and production are on average zero, the pressure 
fluctuations grow diffusively with time. We coin the term 
- diffusive jitter of pressure fluctuations to describe this 
effect. 

• We show that the diffusive jitter of pressure is a nonlocal 
phenomenon where the pressure swings at one location 
depend on the behavior at all other locations. 

• We show that diffusive jitter is spatially inhomogeneous 
and dependent on the spatial distribution of the forecasted 
(stationary) solution. 

The rest of the manuscript is organized as follows. Sections 
[n] and 1^ provide a technical introduction to gas pipeline 
modeling. Section [^provides a brief summary of approaches 
used to solve the steady gas flow problem. Section[V|describes 
a generalization of ID that linearizes the gas flow equations 
around the steady solution that includes the effect of stochastic 
gas consumption. The asymptotic solution of these linearized 
equations describes the diffusive jitter of the pressure fluctua- 



Fig. 1. Schematic representation of the Transco gas transmission network. 


tions. Section [VT| applies the theoretical results to a model of 
the Transco pipeline. Section VII summarizes our main results 
and offers a brief discussion of future work. 


II. Dynamic Gas Flow (DGF) Over a Single Pipe 

Before analyzing a pipeline network, we introduce the gas flow 
equations and notations for a single pipe. Major transmission 
pipelines are typically 16-48 inches in diameter and operate 
at high pressures (e.g. 200 to 1500 psi) and high mass flows 
(millions of cubic feet of gas per day) a, M- Under these 
conditions, the pressure drop and energy loss due to shear 
is modeled by a nearly constant phenomenological friction 
factor /. The resulting gas flow model is a nonlinear partial 
differential equation (PDF) with one spatial dimension x 
(along the pipe axis) and one time dimension im, m, 03: 

3rp-|-3;f(Mp) =0, (1) 

dt{pu)+dx{P'^^)+^xP = -f||^/-pgsina, (2) 

p^pZRT. (3) 

Here, m, p, and p are the spatially-dependent velocity, pressure, 
and density, respectively; Z is the gas compressibility factor; T 
is the temperature, R is the gas constant, and d is the diameter 
of the pipe. 

Eqs. describe mass conservation, momentum balance 

and the ideal gas thermodynamic relation, respectively. The 
first term on the righthand side (rhs) of Eq. Q describes the 
friction losses in the pipe. The second term on the rhs of 
Eq.Q includes the gain or loss of momentum due to gravity 
g when the pipe is tilted by angle a. The frictional losses 
typically dominate the gravity term, which is often dropped. 
Because the flow velocities are usually small compared to the 
sound velocity, the gas inertia term dt{pu) and the advection 
term dx{pu^) are typically small compared to the frictional 
losses and can also be dropped im, ca, ini. For simplicity 
of presentation, we have also assumed that the temperature 
does not change significantly along the pipe. 

Under these assumptions, Eqs. QUID are rewritten in terms 







b) 




Fig. 2. Schematic illustration of a gas network and associated notation, a) 
A schematic illustration of a single edge (i,j) of a network. Nodes at either 
end are indicated by open circles and labeled by their nodal pressure p, and 
pj. Compressors are indicated with filled squares. Mass flow is directed 
from i to j and nodal injections qi and qj contribute to this flow. Nodal 
pressure pj is modified by the compression ratio yielding pijixjj — 0). 
The pressure falls along {i,j} reaching pijixjj —Ljj). If compressor CLj^j is 
not present, then pij(xij =Lij) = pj. b) A schematic of many edges connected 
in a meshed network. Nodes are indexed by / = 0,1,---, where node 0 is 
typically reserved for the swing node - the node where pressure is maintained 
constant throughout the dynamics. Compressors, injections and edge mass 
flows are the same as in a). 


of the pressure p and the mass flux (|) = up: 

cj^d,p + dx<i? = 0, 

P 


2d p 


(4) 

(5) 


where Cs = 'JZRT and p = fZRT are considered constant. The 
solution of Eqs. for t G [0,t] and x G [0,L] requires initial 
and boundary conditions. 


Vx e [ 0 ,L] : (j)( 0 ;x) = (j)o(x), ( 6 ) 

Vf: (j)(f;0) (|)(f;T) (7) 

which are consistent, i.e. (|)o( 0 ) = q'4")(0) and (|)o(T) = 
, in addition fixing the initial pressure at the beginning 
of the pipe, e.g. p(0;0) = po- 


are used to boost pressure. denotes the compression 

ratio of the station adjacent to node / that boosts pressure for 
flow toward node j while OCy^, denotes the compression ratio 
adjacent to node j that boosts pressure for flow toward node i. 
We choose to place compressors at both ends of every edge for 
generality, which also simplifies the notations in the following 
discussion. In reality there is no more than one compressor on 
any particular edge of the graph and a = 1 when there is no 
compressor. Note that 0 C,^j may be larger or smaller than 1, 
allowing the modeling of compression or decompression. If 
only compression is allowed, then > 1. The schematic in 
Fig. [^displays the spatial relationships between nodes, edges, 
and compression ratios. These are expressed mathematically 
as 


Vf G [0 ,t], V(/,;) G E : Pij{t,0) = pi^j{t), (10) 

PijitjLij) — Pj^i{t)j Pi^j ~ Pj^i — Pj^j^ij 

where p, and pi^j are the pressures at node i and pressure after 
compression ratio 0 C,^j . If there is no compressor, then OC,^; = 
1. Under current operating practices, compression ratios do 
not change frequently Thus, we assume that OC,^; does not 
depend on time. 

Eqs. ( |9|10| l are complemented with mass conservation at all 
nodes of the network: 

VfG[0,T], yiG'V: Y. (11) 

When the gas injections q{t) = {qi{t)\i G ‘P’) are given for 
f G [0 ,t], nodal conditions ( [TI) generalizes the single-pipe 
boundary conditions in 0 to a pipe network. Eqs. ® 0 
[T0| [1T| ) constitute a complete set of equations describing the 
Dynamic Gas Flow (DGF) problem when they are supple¬ 
mented with compression ratios, i.e. a= (ai^jj(i,J) G E), 
initial conditions on the flows 


III. Dynamic Gas Flow (DGF) Over a Network 

Next, we generalize the equations for a single pipe to a gas 
network. The network is modeled by a graph £7 = (17, E) with 
vertices 17 and edges E, where the edges are directed (/, J) or 
undirected {i,J} depending on the context. Each vertex i G V 
represents a node with a gas mass injection or consumption 
rate qt. Each edge (i,y) G E is a single pipe with mass flow 
i^ij. The flow along each edge is described by a set of PDFs 
adapted from Eqs. (|4|5|l: 


Vf G [0,t], Vjf,;} G E, Vx G [0;E,;]: 


G %pij{t,x)+d^(ifij{t,x) = 0, 


dxPijitp 




2d 


Pij{t,x) 


= 0 , 


( 8 ) 

(9) 


where pij{t,x) and (|),y(f,x) are the pressure and mass flow, 
respectively, at time f and position x along edge {i,j) of length 
Lij. Here, pij = pjj, (|),j = —(|);„ and L,y =Pji- See Fig. for 
a schematic description of the variables. 

The flow of gas creates a pressure gradient, and compressor 
stations (potentially located at both ends of each edge {i,j}) 


t = 0, G E, Vxij G [Q,Lij]: (|);;(0;x;;) = (j)|“^(x,j),, 

( 12 ) 

and pressure at one slack node, Pi=o(0) = po- 

IV. Optimum Gas Flow Approaches 

Eater in the manuscript, we analyze the DGF problem by 
linearization of the fluctuations about a stationary solution. 
Here, we summarize two approaches to finding this stationary 
solution. We first solve a stationary version of the DGF 
problem, i.e, the Gas Flow (GF) problem, where the steady 
state pressure and flows are expressed in terms of compres¬ 
sion ratios. The time-independent compression ratios are then 
determined via solution of the OGF problem. 

'Compressors are automated to a degree in that they are run in modes 
with one of those modes considered in this manuscript being a constant 
compression ratio. Then pressure fluctuations at the inlet of a compressor 
are amplified at the outlet of that compressor. Other modes of compressor 
operation may lead to different results for the fluctuation amplitude. There is 
relatively fast local control on the compressor to maintain this ratio. However, 
the set points for the ratio is updated very rarely. Adjusting this set point is 
the main effect the pipeline operator has on the system and this adjustment 
is done rai'ely throughout the operating day. 











A. Stationary Gas Flow 

In the GF problem, all input parameters (consump¬ 
tions/injections, compression ratios and the pressure at the 
slack bus) are constant in time. The total injection and con¬ 
sumption is balanced 

(13) 

ieF 

The steady solution of Eq. © is uniform mass flow along 
each pipe in the network, V{/, j} : = const. Substituting 

this result into Eq. (|^ and integrating over space yields the 
algebraic relationship between pressure at position x G [ 0 ;L, j], 
compression, and (constant) flow through the pipe 

V(iJ)€‘E: = pf Ui^j; 

(14) 

The GE problem has a unique solution provided the compres¬ 
sion ratios are known, and the GE solution in ( [T4l l is the basis 
for many approaches to solving the OGE problem. 

B. Optimum Gas Flow 

The solution to the GF problem leaves the time-independent 
compression ratios a unknown. These are chosen by the 
pipeline operators based on a combination of economic and 
operational factors. Here, we describe two approaches for 
selecting the a. The first is a greedy algorithm that approxi¬ 
mates the current pipeline operations in the US. The guiding 
principle is that compressors are activated when the pressure 
prior to the next compressor drops below the acceptable lower 
bound. When activated, a compressor is set to its maximum 
compression ratio a. This algorithm is described in detail 
in 0 . The second approach is based on solutions to the 
optimal gas flow (OGF) problem lfl4l . ifTSll . ifT^ . 171 . Here, 
we summarize a Geometric Programming (GP) approach to 
solving the OGF problem that minimizes the total compressor 
power to move the gas. 

The total power used in pipeline gas compression (assuming 
that the gas is ideal and compression is isentropic) is 

c- 

L (max{a;^^,l}-l), (15) 

where is a constant that depends on the compressor, 
m = (y — 1) /y where y is the gas heat capacity ratio, and 
is the efficiency factor of the compressor. It is important 
to note that fluctuations caused by compressor consumption 
are negligible when compared to gas loads. The term 
denotes the directional mass flow for edge i, j, when the edge 
is oriented from i to j. The OGF formulation assumes that 
the flow through the compressor is from i to j, i.e., > 0 , 

thus the direction of flow must be selected before hand. For 
tree networks, the magnitude and direction of the flows are 
computed exactly apriori and do not depend on the choice of 
compression ratios. For an edge i,j, let ( 7 , and be the two 


(st) 

disjoint graphs obtained by removing {i,j). The flows (^]j 
computed as 




(it) 


L {v) ( 

- - L 

ie Qi ie Qj 


(st) 


( 16 ) 


In networks with loops, flow direction is chosen using heuris¬ 
tics or through the introduction of binary variables na. 
Using the cost function in Eq. the OGF problem is 

formulated as 


mm 

a,/7 


„2 I 


s.t. v(!, ;■)€£: 




= ■ 


Pi 


(17) 


(18) 


ViG'U: 0<p.<pi<Pi, (19) 

V(!, ;■)€£; < ai^j < (20) 

where Eq. ([T8|l is obtained from Eq. ( [T4| . The upper bound in 
Eq. ([g represents engineering limits on pipes and the lower 
bound represents contractual obligations. The upper and lower 
bounds in Eq. ( [20| refer to maximum allowed compression 
and decompression at each compressor. If decompression is 
not allowed, = 1 . 

There are a variety of methods for solving the OGE over trees, 
and in this paper, we use the geometric programming (GP) 
approach described in 0 . The GP approach relaxes the lower 
bound pc,-^ in Eq. (|20jl, i.e. pc,-^, =0. Under this relaxation. 


i-d 


the OGE is transformed into a GP of the form 


mm 

s.t. 


log 


I 44,, 


due" 


Vi e 12 


\( 0 )( 2 E / 

21og(p,.) < Pi < 21og(pi) 

0 < tij < logiUij), 

log < 0 , 


( 21 ) 

( 22 ) 

(23) 

(24) 


V(/,7) G E. 


The transformed variables are related to the original ones via 
the following equations 






(25) 


The OGE in Eqs. ( 2Tp4l l is solved using convex optimization. 
When decompression is not allowed, pti^y=l in Eq. (20i, 


we use a signomial programming (SP) method, which is a 
heuristic version of GP based on solving a sequence of convex 
programs 0 . Here, we use SP to solve the OGF 

V. Diffusive Jitter oe Pressure Fluctuations 
The main contribution of this manuscript builds on the so¬ 
lution of the OGF by introducing a model of stochastic 
gas consumption and the analysis of its effects on pressure 
fluctuations-Vij^Miive jitter. Our approach linearizes the DGF 
equations (Eqs. around a solution to the GE problem 
(augmented with compression ratios from the OGE or the 
greedy compression scenario). 











The linearized model captures the relationship between the 
fluctuating consumption and the fluctuating pressure. Asymp¬ 
totically, the accumulated changes in pressure provide an indi¬ 
cation of how fast the pressure will drift (the jitter) and exceed 
an operating limit in the absence of operator intervention. As 
the DGF solution drifts further from the original GF solution, 
the quality of the linearization degrades. However, we expect 
that the linearized solution to the DGF remains a strong 
relative indicator of how quickly a system will experience 
problems due to stochastic consumption. 

Formally, the stochastic consumption is defined by q{t) = 
q{st) where the components of ^(f) = € 'F') 

are time varying but relatively small in comparison to 
We assume a linearized solution of the DGF problem of 
the form p{t) = -|-5p(f) and (|)(f) = -|- 8 (|)(f), where 

the respective corrections are small, i.e. |5p(f)| ^ and 
|5(|)(f)| <C The linearized versions of Eqs. ([^ |9 10 111 


Vf G y{i,j} G E, Vx G [ 0 ;L, 7 ] : 


c,^dtdpij + dx5(i?ij = 0, 
^x^Pij + ^ I 




2d I ^(st) 


ipfr 


= 0 , 


(27) 


VfG[0,T], y{i,j)e‘E: 

hpi^j = hpiCXi^j^ (28) 

hpij{t,0) = hpi^j{t), hpij{t,Lij) = 5pj^i{t), (29) 
VfG[0,T], ViG-I/: Y. (30) 

We seek asymptotic solutions to the PDF of Eqs. ( |26l [27] 
|28|29|30| l, where asymptotic implies finding solutions for 
time X longer than the correlation time of the fluctuation 
consumption In addition, we seek solutions of Eqs. (26 
27|28|29|30|) that connect the nodal quantities by algebraic 


5pij = aij{t)Zij{x)+bij{t,x), 


(31) 


Substitution of proposed solution ( |3T] ) into Eqs. (26 2^ yields 
an equation for Zij, i.e. 


3^ 


Zij=0, 


(32) 


where Zjj{x) counts x from node i. The integration of Eq. (32i 
over the spatial dependence of the stationary profile 
yields 


Zij{x) — 


2 p!f W ’ 


(33) 


where the normalization constant is chosen to guarantee, 
!QZij{x)dxlL=\. 

We solve for the time-dependent factor «,;■(?) by substituting 
bpij ^ aij{t)Zij{x) into Eq. (26 1 and integrating the result over 


the entire spatial extent of the pipe {/,j} yielding 

aij{t) = mijit’,0) - Hijit\L)). (34) 

(26) In the asymptotic limit where dpij ^ ajj{t)Zij{x) for every 
pipe (graph edge), Eqs. (291 can only be satisfied if the aij{t) 
have the same functional dependence on time, i.e., 

y{i, aij{t) = a{t)cij, (35) 

where Cij = cji is an edge specific constant. 

To compute the global time-dependent factor a{t) we sum the 
mass conservation equation over all the nodes of the graph 

Y^i= E (36) 

integrate over time and define 

m= fdt'YW), 


(37) 


relationships thereby eliminating the complexity of the orginal 
PDE. 

The solution approach is an extension of the work in |[8l . 
Following mi, we solve Eqs. ( |26| jZTj ) for each pipe using 
a proposed solution of the form 


and finally sum Eq. (351 overall edges; 

L{ij}e'ECij 

Therefore, Vf, y{i,j}G‘E, x G [0,Lij] : 


hpij{t,x) 




(38) 


(39) 


where the aij{t) depend on time. In O, it was argued that 
the aij{t)Zij[x) term represents the asymptotic contribution 
to the gas pressure fluctuations that grows in time. In con¬ 
trast, bij{t,x) represents smaller contributions to the pressure 
fluctuations that do not grow in time. Here, we focus on the 
contribution from the a,j(f)Z, 7 (x) term which is asymptotically 
dominant at long timesH 

^Bounds for the second term are derived and solved using inhomogeneous 
linear equations for bjj. 


T.{iJ}e‘E Cij 

The unknown edge constants ctj are derived by substituting 
Eqs. ( |39| l into Eqs. ( [28] [35] l yielding 

V/, V;,k s.t. G E : ( 40 ) 

Eqs. ( |39||40l[33] ) express the complete asymptotic (zero mode) 
solution of the DGF problem. 

Finally, we make several observations to connect the solution 


for the pressure fluctuations in Eqs. (39 40 33 1 to a prob¬ 
ability distribution over the pressure fluctuations. First, the 
random gas load fluctuations ^, (f) are zero-mean, temporarily 
homogeneous, and relatively short correlated in both time (the 
correlation time is less than x) and space (the correlation length 



































is less than the spatial extent of the network). Second, the 
fluctuations of bpij in Eq. (39i are given by a time-integral 


and spatial-sum of the fluctuations. According to the Large 
Deviation theory, these observations imply that the pressure 
fluctuations form a Gaussian random process which jitters 
diffusiveiy in time. Specifically, the Probability Distribution 
Function (PDF) of dpij{t,x) is 


T{bpij{t,x) = 5 ) {2ntDij{x)) 


-1/2 


exp I 


2tDij (jr) 


Dij = 


f clcjjZjflx) 
\T.{k,l}e'ECkl 



E Ut') 

neV 


E Ut') 



*N. 


estimate of the gas consumption fluctuations assumes that 
the fluctuations at neighboring nodes are uncorrelated. If 
these neighboring nodes are gas-fired turbine generators that 
are both being used to balance renewable fluctuations, the 
assumption of independence may lead to an underestimation 


in Eq. (431. 


,(41) 


(42) 


where the correlation function on the right-hand-side does not 
depend on f due to assumption of the statistical homogeneity 
of t 

VI. Numerical Experiments 


Inspection of Eq. (42 1 shows that the variance of the pressure 


fluctuations as a function of position in the network is related 
to the coefficients Diflx), referred to collectively as D. Higher 
values of D correspond to larger pressure fluctuations and 
higher likelihood of the pressure violating an engineering or 
contractual limit. By analogy with related physical processes, 
the coefficients D are similar to a diffusion coefficient, and we 
refer them this way in the remainder of the manuscript. The 
origins of D are primarily twofold. Once the gas consumptions 
and injections are fixed, the spatial dependence of D arises 
from the particular stationary solution of pressures, flows, and 
compression ratios through the Zijix). The magnitude of D is 
also related to the average global strength of the consumption 
fluctuations 

We apply the results described above to the Transco pipeline 
shown schematically in Fig. [T] We use data for the total con¬ 
sumption at each node over a 24-hour period from December 
29, 2012 to fix the forecasted consumption for the stationary 
GF solution. These data represent relatively stressed operations 
for the Transco pipeline. The Transco pipeline has a small 
number of loops, which we partition to create a tree topology 
0 that is very nearly linear but with a few small branches. 
We resolve these branches in the solution of the GF (or OGF) 
problem, however, when analyzing the pressure fluctuations, 
we aggregate these short branches to nodal consumptions and 
only analyze the fluctuations as a function of distance along 
the mainline. 

The Transco operational data does not include information 
on the deviations of the gas flows from their average or 
scheduled values. Instead, we estimate the global mean-square 
consumption fluctuations as 


For presentation purposes, it is convenient to find a suitable 
normalization for D. Motivated by Eqs.( 41|431 l, we normalize 
D by 

^o«(y) Ao 

where po = 800 psi « 5.5 * 10® Pa is the upper bound on 
allowed pressure in the pipes and to = 15 min « lO^s is a 
representative time period where we expect the developed 
theory to work well. 

We consider the base case of December 29th, 2012 and 
several modifications of this base case to investigate the 
effects of changing operations. Fig. [^displays D as a function 
of location along the mainline for two different stationary 
solutions for the base case—the OGF solution described in 
Section IV and the greedy algorithm from 0. For a charac¬ 
teristic time of 15 min « lO^i, D/Dg = 1 corresponds to a 
variance in pressure fluctuations of (266 psi)^. Pressures in 
the Transco Pipeline range between 500 psi and 800 psi, so 
the pressure fluctuation standard deviation is 33 — 53% of the 
pressures in the pipeline for D/Dg = 1. The same characteristic 
time and D/Dg = 0.1 yields a variance of (84 psi)^ which 
gives pressure fluctuation standard deviations of 10— 16% of 
pipeline pressures. Since the pressure variance grows linearly 
in time, often over several 15 min intervals, these fluctuations 
can quickly grow to exceed pressure bounds without proper 
intervention. As the plots show, most pressure fluctuations are 
above DjDg =0.1 throughout the pipeline, and therefore the 
fluctuations are of concern in any regions of pipeline where the 
pressure is near its upper or lower bound. The two solutions 
display similarities. Both show a build up of D from milepost 
800 nearer to the Gulf of Mexico, a peak at milepost 1771 near 
New York and New Jersey, and a decay to a smaller value at 
milepost 2000 near the injection point for the Marcellus Shale 
in Pennsylvania. 

The spatial variation of D is due to Ziflx), and the general 


(43) 


Here, (|)o « 20 kg/s is a typical average consumption for a 
node in the Transco pipeline, and V « 70 is the number 
of consumption nodes, e.g. city-gates or power plants. This 


shape of D can be understood by revisiting Eq. 32 The 
form of this equation suggest exponential growth or decay 
of Zij(x) depending on the orientation of (|),y. The flow from 
the Gulf to the New York/New Jersey area is unidirectional 
creating the growth of D observed in Fig. However, the 
large loads in the New York/New Jersey area combined with 
the offsetting injections from from the Marcellus Shale creates 
a flow reversal and an exponential decay of Zij (and therefore 
of D). The peak in D is connected to the point of flow reversal. 
The solution in Fig. displays more structure than simple 
exponential growth and decay for several reasons. First, the 
mass flow rates (|),y depend on location. However, perhaps 
more important are the discontinuities in D. These occur at 
compressor stations and are due to the discontinuities in p^^^^ 












at these locations. 

The global behavior of the OGF solution and the greedy 
algorithm Fig. is similar because the differences in the 
compression ratios in the stationary solution does not affect the 
mass flow rates. However, it does affect the spatial dependence 
of pressure which can lead to the substantial local differences 
observed in Fig. For example, between mileposts 1400 and 
1700, D is much lower for the greedy algorithm compared 
to the OGF. However, the deployment of a compressor near 
milepost 1700 in the greedy algorithm leads to a large jump 
in D, a larger peak in pressure fluctuations, and a greater 
chance for violation of an engineering or contractual pressure 
limit. For this one example, this difference would seem to 
suggest that the OGF solution is less susceptible to pressure 
fluctuations. However, we note that the expected magnitude of 
the pressure fluctuations is not taken into account in either the 
greedy algorithm or the OGF. 

What these results do suggest is that the deployment of 
compressors in the stationary solution can have a signihcant 
impact on the expected pressure fluctuations, and that it is 
possible to formulate a a compressor dispatch optimization 
that balances the risk of such fluctuations against other desired 
operational properties, e.g. cost. The simple algebraic form of 
the probability of such large fluctuations in Eqs. ( 41|42 i are 
convenient for incorporation into such formulations. 

To determine the effect of overall consumption and injection 
on D, we uniformly scaled the base case consumptions and 
injections by a constant factor—a scaling that preserves the 
balance of consumptions and injections required for the exis¬ 
tence of a stationary solution. Figure displays the results for 
D computed using the OGF stationary solutions. The stationary 
solutions show small local differences in D caused by the 
deployment of compression. However, the major impact stems 
from the increase (or decrease) in mass flows. The uniform 
scaling does not affect the location of the flow reversal, so 
the peak in D appears at the same place. However, the larger 
(smaller) flows lead to faster (slower) growth rates for Z,y (see 
Eq. [3^ and an overall higher (lower) peak in D. 

In recent years the Marcellus Shale has become a large 
supplier of gas, and its injection capability is expected to 
increase El. To model the effect of this expansion we scaled 
all injections from the Marcellus Shale by a constant factor and 
removed a corresponding amount of gas from the injections 
at the Gulf to preserve the global balance of consumption 
and injection. Eig. displays the results for D along the 
Transco pipeline for the OGE solution. Although the injection 
from the Marcellus is increased, the major gas load centers 
in New York and New Jersey keep the flow reversal point, 
and therefore the peak in D, pinned at more or less the same 
location. Larger Marcellus injections show slightly lower peak 
amplitudes of D indicating that, as gas injections are shifted 
from the Gulf to the Marcellus Shale, the reliability of pipeline 
operations is improved. We conjecture that this change is due 
to moving the source of gas injections closer to the major load 
centers. However, we again note that the OGE methods used 
to And the stationary solution do not account for the expected 



Fig. 3. Diffusion coefficient as a function of distance along the Transco 
mainline with stationary solutions given by OGF and the greedy algorithm. 
Both show a peak at milepost 1771, but the magnitude of this peak is much 
higher for the greedy algorithm than the OGF, indicating larger pressure 
fluctuations for the greedy algorithm. 


pressure fluctuations and their inclusion will likely lead to a 
modification of these results. 

Although scaling overall consumption and Marcellus supply 
are directly relevant to operators and planners, they do not 
exhibit a shift in the location of the peak in D or the appearance 
of multiple local maxima. To study these possibilities, we 
next imposed some less realistic changes. In particular, the 
loads in New York are shifted to points closer to the Gulf and 
Marcellus Shale, but the New Jersey loads were left unaffected. 
Eig. 1^ displays the impact of this shift on D computed using 
stationary solutions from the OGE and the greedy algorithm. 
Since the large load in New Jersey remains, milepost 1771 is 
still a position of flow reversal and a, now minor, maximum 
of D. However, the redistribution of load leads to a new 
global maximum near milepost 1319—a large load in North 
Carolina. Although the location of maximum fluctuations has 
been relocated, the maximum of D is much reduced by moving 
the loads closer to the gas injections. 

The previous example showed the appearance of a new global 
maximum, as well as several small local maxima, but the orig¬ 
inal local maximum remained. To remove it, we shifted New 
Jersey’s load to points closer to the Gulf and the Marcellus 
Shale, as shown in Eig. This successfully removes the local 
maximum at milepost 1771 while leaving the global maximum 
at milepost 1319. In this case the jitter (diffusion coefficient) 
of the greedy algorithm and OGE are comparable, with OGE 
jitter greater before milepost 1319 and greedy greater after 
milepost 1319. 

VII. Conclusions and Euture Work 

We have focused the analysis on the coupling between natural 
gas networks and electric networks at the time scale of intra¬ 
day natural gas markets. The coupling at this time scale is 
expected to become tighter because of several factors: the 
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Fig. 4. Diffusion coefficient as a function of distance along the Transco 
mainline with stationary solutions given by the OGF with global consumption 
and injection scaled by a uniform factor. All show a peak at milepost 1771, 
but higher scaling factors have higher magnitudes at their peaks, indicating 
larger pressure fluctuations for larger system loads. 
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Fig. 6. Diffusion coefficient as a function of distance along the Transco 
mainline with load redistributed from the large load in New York to the Gulf 
and Marcellus Shale, leaving the large load in New Jersey unaltered. This 
causes the appearance of a new global maximum at milepost 1319 which is 
the location of a large load in North Carolina. Since the New Jersey load was 
not redistributed, a local maximum remains at milepost 1771. 
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Fig. 5. Diffusion coefficient as a function of distance along the Transco 
mainline with stationary solutions given by the OGF with Mai'cellus Shale 
injections scaled by a factor and the corresponding amount of injections 
removed from the Gulf. All show a peak at milepost 1771, but higher scaling 
factors have slightly lower magnitude at the peak, indicating smaller pressure 
fluctations. Higher scaling factors also have much lower magnitudes in the 
Marcellus Shale, indicating smaller pressure fluctuations when injections are 
shifted from the Gulf to the Marcellus Shale. 

retirement of coal and fuel oil-fired generation in favor of 
natural gas-fired generation because of environmental concerns 
and the increased availability and low cost of natural gas and 
the ability of gas-fired generation to respond quickly to the 
variability of renewable generation. Although gas pipelines 
have the ability store gas in the form of increased pressure 
in the piepline, i.e. linepack, this storage is limited. In the 
future, linepack will be increasingly exercised as more gas- 
fired generation is used to balance increasing amounts of wind 
generation. Larger swings in gas pipeline pressure (linepack) 
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Fig. 7. Diffusion coefficient as a function of distance along the Transco 
mainline with load redistributed from the large loads in New York and New 
Jersey closer to the Gulf and Marcellus Shale. The global maximum at 
milepost 1319 remains, but the local maximum at milepost 1771 dissappears 
since the large load has been removed from this area. 

affect the ability of the pipeline to deliver gas to the generators, 
creating reliability implications that cascaded across these two 
infrastructures. 

In this initial work, we have assessed the impact of fluctuating 
consumption by gas-generators on pipeline pressure. We start 
by splitting the gas flow equations for pipelines into two parts. 
The first is a stationary part that is time-independent and 
reflects the gas flows scheduled by the gas markets and gas 
compressor deployment determined by the pipeline operator. 
The second is a representation of the fluctuations around the 
scheduled flows created by linearizing the gas flow equations 
about the scheduled flows and compressor operations. From 
this linearized model, we can predict the probability that a 
set of stochastic gas loads will cause the pipeline pressure to 
violate an engineering or contractual pressure limit and create 









































a reliability concern for the pipeline operator or the electrical 
grid operator. By making assumptions about the nature of the 
gas consumption fluctuations, this probability can be expressed 
in an algebraic form that is convenient for integration into 
a gas flow/gas compressor optimization problem where the 
probability can be a constraint or part of the objective to limit 
the likelihood of a pipeline reliability issue. 

We applied the theoretic results to a realistic model based on 
the Transco pipeline. Our computational experiments with the 
Transco model revealed the following interesting observations. 
First, the probability of large pressure fluctuations is highest 
at locations in the pipeline where the gas flow experiences a 
reversal (in and around the New York/New Jersey area for the 
Transco pipeline). Second, increasing the stress on the pipeline 
by increasing gas flow rates leads to higher probabilities of 
large pressure fluctuations. Third, rearranging pipeline flows, 
e.g. by increasing purchases from the Marcellus Shale at the 
expense of gas from the Gulf, can decrease the probability of 
large pressure fluctuations by moving the gas sources closer 
to the gas loads. 

The results of this paper suggest a number of interesting 
directions for future research. 

• The linearization and asymptotic assumptions described 
here need to be validated against direct dynamic (tran¬ 
sient) simulations of gas flows in variety of situations. 
Most existing work on such validations ca, ca, im, 
m, im, Ea uses single pipe models. The challenge 
is to develop fast computational algorithms for transient 
problems with mixed (initial and boundary) conditions 
over large and loopy gas networks. 

. Our dynamic method applies to gas networks with loops, 
back flows, bi-directional compression and other compli¬ 
cations. We will extend the experimental study to other 
current and planned networks in the U.S. and Europe. 

• Extending the probabilistic risk framework to the compli¬ 
cations mentioned above requires extending the methods 
of Cl to create efficient optimization algorithms for gas 
networks with loops. 

• Compressor positions are assumed fixed by the OGE 
and greedy solution methods presented here. However, 
compressor position has a significant effect on pressure 
fluctuations in its vicinity. A future direction will be 
to formulate a compressor dispatch scheme and use it 
to analyze the effect of varying compressor position on 
pressure fluctuations near the compressors. 

• Incorporation of the probabilistic risk measures into OGE 
formulations to directly account for this risk. A promising 
direction is the chance constrained methodology devel¬ 
oped in E^ . 

• This work suggests a new mathematical, statistical and 
computational foundations necessary to address the com¬ 
prehensive strategic problems of re-organizing the exist¬ 
ing system of energy trading (in the U.S. and elsewhere). 
Such a reorganization is required to reduce inefficiencies 
in how power and gas markets interact m, El, EH- 

• It is not realistic to expect (at least not in US) that 


gas and power markets will merge in the near future. 
However, it is important to account for effects of mutual 
dependencies. In particular, incorporating effects of gas 
pressure fluctuations and uncertainty into planning and 
operations of power systems with significant penetrations 
of renewables and with gas turbines involved in balancing 
the renewable fluctuations is a very promising future 
direction for research. On the other hand it is as important 
to account for the effect of ramps in gas consumptions at 
generators on the gas flow optimization. 
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